/**************************************************************/
/*  This file is licensed from, and is a trade secret of:     */
/*                                                            */
/*    Mercury Computer Systems, Incorporated                  */
/*    199 Riverneck Road                                      */
/*    Chelmsford, Massachusetts 01824-2820                    */
/*    United States of America                                */
/*    Telephone + 1 978 256-1300                              */
/*                                                            */
/*  Refer to your Software License Agreement for provisions   */
/*  on use, duplication, or disclosure.                       */
/*  US. Government users and prime contractors are also       */
/*  governed by the following Data Rights Notice:             */
/*  Use, duplication or disclosure by the US. Government is   */
/*  subject to the provisions set forth in FAR 52.227-14 for  */
/*  restricted rights or other agency data rights provisions  */
/*  on commercial software and documentation, as applicable.  */
/*                                                            */
/*  This file: (a) was developed at private expense, and no   */
/*  part was developed with government funds, (b) is a        */
/*  trade secret of Mercury Computer Systems, Inc. for the    */
/*  purposes of the Freedom of Information Act, (c) is        */
/*  "commercial computer software" subject to limited         */
/*  utilization as provided in the above noted Software       */
/*  License Agreement and (d) is in all respects data         */
/*  belonging exclusively to Mercury Computer Systems, Inc..  */
/*                                                            */
/**************************************************************/

/******************************************************************/
/* Copyright (c) 1999-2003, Mercury Computer Systems, Inc.,       */
/* Chelmsford MA.,and all third party embedded software sources.  */
/* All rights reserved under the Copyright laws of US. and        */
/* international treaties.                                        */
/*                                                                */
/* Patent Pending                                                 */
/******************************************************************/


//$COPYRIGHT_END

#include <string>
#include <cstdlib>
#include <cstdio>
#include <cassert>
#include "RsHop.h"

using namespace std;

extern "C" {
  RSHOP_EXPORT
  POA_CF::Resource *
  CreateRsHopDecoder (const char * identifier)
  {
    RsHopDecoder * decoder;

    if (RsHopDecoder::Create (decoder, identifier) != success) {
      return 0;
    }

    return (POA_CF::Resource *) decoder;
  }
}

static const unsigned long max_output_buffer_size = 168;

status
RsHopDecoder::Create (RsHopDecoder *& handle,
                      const char * identifier)
{
  RsHopDecoder * Temp = new RsHopDecoder (identifier);

  if (!Temp) {
    return fail;
  }

  if (!Temp->Init ()) {
    handle = 0;
    return fail;
  }

  handle = Temp;
  return success;
}

RsHopDecoder::RsHopDecoder (const char * identifier)
  : CF__Resource_impl (identifier),
    m_threadrunning (false),
    m_bStop (false)
{
  unsigned int input_buffer_size;

  // Initialize buffers.
  input_buffer_size = 176; // Able to handle 4 buffers @ 16 bytes from FM3TR receiver.
  m_input = new char[input_buffer_size];
  m_output = new char[max_output_buffer_size];

  // Define ports.
  m_rshopin = new RsHopInPort_impl (input_buffer_size);
  m_rshopout = new BitPacketizerOutPort_impl (m_output, max_output_buffer_size);
}

RsHopDecoder::~RsHopDecoder ()
{
  delete m_input;
  delete m_output;
  delete m_rshopin;
  delete m_rshopout;
  delete m_gf;
  delete m_decMsg;
}

status
RsHopDecoder::Init ()
{
#ifdef MCS_VERBOSE
  std::cout << "RsHopDecoder: Initializing." << std::endl;
#endif

  unsigned int hop_prim_poly[6] = {1, 0, 1, 0, 0, 1};

  // Initialize the parameters.
  m_n = 16;
  m_k = 14;
  m_N = 31;
  m_K = 29;
  m_m = 5;
  m_b = 1;

  m_gf = new galois::GaloisField (m_m, hop_prim_poly);
  m_decMsg = new galois::GaloisFieldPolynomial (m_gf, m_N-1);

  if ((m_n-m_k)%2==0) {
    m_t = (m_n - m_k)/2;
  }
  else {
    m_t = (m_n - m_k - 1)/2;
  }

  return success;
}

void
RsHopDecoder::RunWorker (void * param)
{
  RsHopDecoder * component =
    (RsHopDecoder *) param;

#ifdef MCS_VERBOSE
  std::cout << "RsHopDecoder: Component is running." << std::endl;
#endif

  component->RunLoop ();

#ifdef MCS_VERBOSE
  std::cout << "RsHopDecoder: Component is stopping." << std::endl;
#endif
}

void
RsHopDecoder::start ()
  throw (CORBA::SystemException, CF::Resource::StartError)
{
  if (m_threadrunning) {
    return;
  }

  m_threadrunning = true;
  m_thread.start (RunWorker, (void *) this);
  CF__Resource_impl::start ();
}

void
RsHopDecoder::stop ()
  throw (CORBA::SystemException, CF::Resource::StopError)
{
  m_bStop = true;
}

CORBA::Object_ptr
RsHopDecoder::getPort (const CORBA::Char * name)
  throw (CORBA::SystemException, CF::PortSupplier::UnknownPort)
{
  string temp (name);

  if (temp == "RsHopDecoderIn") {
    return m_rshopin->getPort ();
  }
  else if (temp == "RsHopDecoderOut") {
    return m_rshopout->getPort ();
  }
  else {
    throw CF::PortSupplier::UnknownPort ();
  }
}

status
RsHopDecoder::RunLoop ()
{
  m_bStop = false;
  m_output_count = 0;

  while (!m_bStop) {
    Run ();
  }

  return success;
}

status
RsHopDecoder::Run ()
{
  unsigned int i, num_read=0;
  unsigned int m_input_count = m_rshopin->GetData (&m_input);

  if (m_input_count == 0) {
    return fail;
  }
  else if (m_input_count < m_n) {
    // Not enough data in buffer; release and wait for more.
    m_rshopin->ReleaseDataBuffer (0);
    return fail;
  }

#ifdef MCS_VERBOSE
  std::cout << "RsHopDecoder: Received " << m_input_count
	    << " samples." << std::endl;
#endif

  while (m_input_count >= m_n) {
    // Initialize the gf polynomial with the incoming buffer.
    for (i=0; i < m_n; i++) {
      (*m_decMsg)[i] = m_input[i];
    }

    /* Optionally, add errors to the signal. */

    galois::GaloisFieldPolynomial hopErrors (m_gf, m_N-1);

    /* Calculate syndromes. */

    // initialize the syndrome vector with zeros
    galois::GaloisFieldPolynomial S_hop (m_gf, m_n - m_k - 1);
    bool syndromeZeroFlag = false;

    for (i=0; i< m_n - m_k; i++) {
      S_hop[i] = (*m_decMsg) (m_gf->alpha (m_b+i));
      syndromeZeroFlag = syndromeZeroFlag || (S_hop[i].poly() == 0);
    }

    /*
     * Calculate error locator L(x) and magnitude W(x) polynomials.
     *
     * Do this step only if Syndrome vector is non-zero:
     * If first element of syndrome vector (polynomial form) is zero,
     * then skip finding error locator and magnitude polynomials.
     */

    if (syndromeZeroFlag) {
      for (i=0; i<m_k; i++) {
	m_output[m_output_count + i] = m_input[i + m_n - m_k];
      }

      m_input_count -= m_n;
      num_read += m_n;
      m_input += m_n;
      m_output_count += m_k;
      continue;
    }

    /*
     * There is a bug with the galois library, such that is one of the
     * two polynomials is of degree 1, then the polynomial division does
     * not return anything.  Therefore, both the numerator and denominator
     * polynominals are multiplied by 'x' for the time being.  This has
     * no effect on the remainder and the quotient, and is basically
     * adding a 0 to both the denominator and the numerator.
     */

    galois::GaloisFieldPolynomial A_1 (m_gf, 2*m_t);
    galois::GaloisFieldPolynomial A_2 (m_gf, 2*m_t - 1);

    galois::GaloisFieldPolynomial D_1 (m_gf, 2*m_t);
    galois::GaloisFieldPolynomial D_2 (m_gf, 2*m_t);

    // Dividend set to x^2t initially.
    // A_1: {0,0,1}

    for (i=0; i<2*m_t; i++) {
      A_1[i] = 0;
    }

    A_1[2*(m_t)]  = 1;

    // A_2 = Syndromes
    for (i=0;i<2*m_t;i++) {
      A_2[i]  = S_hop[i];
    }

    // Initialize D polynomials.
    // D_1 is already set to all zeros.
    D_2[0] = 1;

    unsigned int degree = 2*m_t;

    galois::GaloisFieldPolynomial Q_temp, R_temp, D_temp;

    while (degree >= m_t) {
      /*
       * Add one degree to A_1 and A_2 in order to work around the
       * algorithm problem.
       */

      A_1 <<= 1;
      A_2 <<= 1;

      /*
       * The line below is the polynomial division and the most
       * cycle-intensive part of the algorithm.
       */

      Q_temp = A_1 / A_2;
      R_temp = A_1 % A_2;

      // Remove the degree added to work around the algorithm problem.
      R_temp >>= 1;

      // Update A_1 and A_2
      A_1 = A_2;
      A_2 = R_temp;

      // Get the polynomial degree of R_temp, to check for the while loop.
      degree = R_temp.deg();

      // Update the D polynomial.
      D_temp = (Q_temp*D_2) + D_1;

      // Update D_1 and D_2
      D_1 = D_2;
      D_2 = D_temp;
    }

    D_2.simplify();

    if (D_2[0] != 0) {
      if (D_2.deg()<=m_t+1 && D_2.deg()>=1 && A_2.deg() <= m_t) {
	// W polynomial from A_2
	A_2 = A_2 / D_2[0];

	// We've got the L polynomial from D_2.
	D_2 = D_2 / D_2[0];
      }
      else {
	// We have a failure in decoding the signal.
	D_2 = 0;
	A_2 = 0;
      }
    }
    else {
      if (D_2.deg()<=m_t+1 && D_2.deg()>=1 && A_2.deg() <= m_t) {
	// Do nothing.
      }
      else {
	// We have a failure in decoding the signal.
	D_2 = 0;
	A_2 = 0;
      }
    }

    /* Calculate the error locations. */
    unsigned int * errorLocations = new unsigned int [m_N];
    unsigned int numLocations = 0;

    /*
     * Solve the error locator polynomial (by chien search).  To find the
     * error locator polynomial roots, basically try all of the alpha
     * values one by one. This is a trivial search since one needs to
     * try `n` number of alpha values (from alpha^0 ... alpha^(n-1)).
     *
     * L is evaluated at [alpha^(-1) ... alpha^(-14) alpha^(-15)].
     * Zeros give the roots of the L (error locations).
     *
     * For N>n (Shortened codes), if an error location is greater than
     * n, it should be interpreted as a decoding error.
     */

    for (i=1;i<m_N+1;i++) {
      if (D_2 (m_gf->alpha((m_N)-i)) ==0 ) {
	errorLocations[numLocations] = i;
	numLocations++;
      }
    }

    /*
     * Solve the error magnitudes, and correct the errors (using
     * Euclidean Algorithm).  For euclidean algorithm formulation,
     * check out:
     *   http://en.wikipedia.org/wiki/Extended_Euclidean_algorithm
     *   http://www.dsprelated.com/showmessage/22571/1.php
     */

    // First, obtain L'(X_j^(-1)) = L_1 + L_3*X_j^(-2) + L_5*X_j^(-4) + ...
    galois::GaloisFieldPolynomial L_prime (m_gf, D_2.deg());

    // Remove the constant term.
    for (i=0; i<D_2.deg() ; i++) {
      L_prime[i] = D_2[i+1];
    }

    // Set even terms to zero.
    for (i=1; i<L_prime.deg()+1; i+=2) {
      L_prime[i] = 0;
    }

    L_prime.simplify();

    galois::GaloisFieldElement errorVal (m_gf,0);
    galois::GaloisFieldPolynomial hopDecoderOut = *m_decMsg;

    for (i=0; i<numLocations; i++) {
      if (L_prime (m_gf->alpha (m_N-errorLocations[i])) != 0) {
	errorVal  = m_gf->alpha ((errorLocations[i]*(1-m_b)) % m_N);
	errorVal *= A_2 (m_gf->alpha (m_N - errorLocations[i]));
	errorVal /= L_prime (m_gf->alpha (m_N - errorLocations[i]));
	hopDecoderOut[errorLocations[i] % m_N] += errorVal;
      }
      else {
	// Decoding of the codeword failed.
      }
    }

    for (i=0; i < m_k; i++) {
      m_output[m_output_count + i] = (char) hopDecoderOut[i + m_n - m_k].poly();
    }

    m_input_count -= m_n;
    num_read += m_n;
    m_input += m_n;
    m_output_count += m_k;
    delete errorLocations;
  } // while

  if (m_output_count >= 140) {
#ifdef MCS_VERBOSE
    std::cout << "RsHopDecode: Sending " << m_output_count
	      << " samples." << std::endl;
#endif
    assert (m_output_count <= max_output_buffer_size);
    m_rshopout->pushData (m_output, m_output_count);
    m_output_count = 0;
  }

  // Release the input buffer
  m_rshopin->ReleaseDataBuffer (num_read);
  return success;
}
